\chapter{Signal Processing Tools}
\label{cha:signalprocessing}

      This chapter presents the various signal processing tools
  provided by the library.  

\section{1D Transforms}
\label{sec:1dtransforms}

    All transforms in \libit share the same framework. They are
  objects taking a generic vector (Vec) as input and output a
  transformed generic vector, which may be of different type or length
  (e.g. for redundant transforms). Declaring new transforms is done by
  inheriting from the it\_transform\_t object and defining the transform
  and inverse transform methods. For more details on objects and
  object oriented programming in C with libit, see
  Chapter~\ref{cha:objects}. Here is the set of transforms provided by
  libit (which is currently quite limited). 

\section{Discrete Fourier Transform}
\label{sec:dft}

    The Fourier transform of a signal gives a frequency representation
  of the energy and relative phase present in this signal. The
  discrete transform \emph{X} of a signal
  \emph{x} is defined in complex domain as: 

$$
X[n] = \sum\limits_{k=0}^{N-1} x[k] e^{-2i\pi \frac{kn}{N}}.
$$

    This transform is invertible and the original signal can be
    recovered by using the inverse transform:

$$
x[k] = \frac{1}{N} \sum\limits_{n=0}^{N-1} X[n] e^{2i\pi \frac{kn}{N}}.
$$


    In libit, the discrete Fourier transform of a vector is obtained
  by calling \function{it\_dft()}, which takes a
  complex vector (cvec) as input and returns a complex vector
  representing the amplitude and phase for each frequency. The inverse
  transform is obtained by using \function{it\_idft()}.

    Currently the transform is implemented as in the math formula,
  which is inefficient ($O(N^2)$ vs $O(N \log N)$ for the FFT
  algorithm). However it has the advantage being able to handle
  vectors of any length. For vectors where the length can be
  decomposed in a power of two times an odd number, the current
  implementation of the DFT could be optimized much. Also, there is no
  optimized function for the special case of real samples yet; Thus,
  real-valued vectors (vec) must be converted to complex vectors using
  \function{vec\_to\_cvec()} before calling
  \function{it\_dft()} (which will lead to a
  symmetric conjugate frequency representation). 

\begin{program}{Discrete Fourier transform example}{pro:dft}
/* analyse the vector using the discrete Fourier transform */
cvec vt = it_dft(v);

/* recompose the vector by inverse transform */
cvec vr = it_idft(vt);
\end{program}

\section{Discrete Wavelet Transform}
\label{sec:discretewaveletransform}

    The dicrete wavelet transform (DWT) analyses a signal into a given
  number of subbands (scales) while maintaining a time-frequency
  localization compromise. High frequencies and impulses are well
  localized in space, while low frequencies are precise. The DWT
  provided by libit is a dyadic, reversible transform implemented
  using the lifting scheme. Currently there is no factorization
  algorithm in the library to extract the lifting steps from the
  definition of the FIR filters. Therefore these steps must be given
  explicitly. However, the most commonly used 5/3 and 9/7 biorthogonal
  wavelet lifting steps are provided.  

    To apply a wavelet transform on a vector, use the
  \function{it\_dwt()} function. It takes a real
  vector as input and output a real vector of the same length
  containing all the concatenated subbands (from low-frequency to
  high-frequency). Additional parameters include the number of levels
  and the lifting steps used for the decomposition. The output vector
  can be split in an array of subbands using
  \function{it\_wavelet\_split()}, and merged back using 
  \function{it\_wavelet\_merge()}. The
  inverse transform is provided through 
  \function{it\_idwt()}. Here is an example on how to
  analyse a signal using the 9/7 biorthogonal wavelet. 

\begin{program}{Discrete wavelet transform example}{pro:discretewavelettransform}
/* analyse the vector using a 5-level 9/7 wavelet decomposition */
vec vt = it_dwt(v, it_wavelet_lifting_97, 5);

/* clear the high subband */
vec_set_between(vt, (vec_length(v) + 1) / 2, end, 0);

/* recompose the vector by inverse transform */
vec vr = it_idwt(vt, it_wavelet_lifting_97, 5);
\end{program}

\section{2D Transforms}
\label{sec:2dtransforms}

    Similar to the 1D case, 2D transforms share the same
  framework. They are objects taking a generic matrix (Mat) as input
  and output a transformed generic matrix, which may be of different
  type, width or height as the original matrix (e.g. for redundant
  transforms). Declaring new 2D transforms is done by inheriting from
  the it\_transform2D\_t object and defining the transform and inverse
  transform methods. For more details on objects and object oriented
  programming in C with libit, see Chapter~\ref{cha:objects}. Here is
  the set of 2D transforms provided by libit (which is currently quite
  limited too). 

\section{Generic Separable 2D Transform}
\label{sec:genericseparable2dtransform}

    A generic 2D separable transform is provided. This transform takes
  a 1D transform and applies it on the row, then the columns of the
  input matrix. Here is an example on how to perform a 2D Fourier
  transform by combining the 1D Fourier transform with the separable
  transform.  

\begin{program}{Separable transform example}{pro:separabletransform}
/* create a Fourier transform */
it_fourier_t fourier = it_fourier_new();
/* create a separable transform out of it */
it_separable2D_t fourier2D = it_separable2D_new(fourier);

/* apply the 2D discrete Fourier transform on a matrix */
cmat mt = (cmat) it_transform2D(fourier2D, (Mat) m);

/* inverse the 2D discrete Fourier transform */
cmat mr = (cmat) it_itransform2D(fourier2D, (Mat) mt);
\end{program}

\section{Separable 2D Wavelets}{pro:separable2dwavelets}

     Bidimensional separable wavelet decomposition is provided using
the lifting implementation. The lifting procedure guarantees perfect
reconstruction and provides an efficient way of filtering. Currently
only the 9/7 and 5/3 wavelets are provided.  

     A new 2D wavelet transform is created from the lifting steps, the
width and height of the image to transform and the number of
decomposition levels. As a particular kind of 2D transform, transform
and inverse transform are provided using the it\_transform2D() and
it\_itransform2D() calls. Here is an example using the 9/7 wavelet:

\begin{program}{2D wavelet transform example}{pro:2dwavelet}
  /* read your favorite version of lena */
  mat m = mat_pgm_read("lena.pgm");

  /* create a new 9/7 wavelet object for a decomposition on 5 levels */
  it_wavelet2D_t *wavelet2D = it_wavelet2D_new(it_wavelet_lifting_97, 5);

  /* transform the image */
  mat mt = it_transform2D(wavelet2D, m);

  /* reconstruct the image */
  mat mr = it_itransform2D(wavelet2D, mt);
\end{program}